clear all 
use "data/both.dta", clear

set seed 2025
* ANOVA results used in note of figure 1. 
oneway roll1 treatment if experiment==1

* Create the results frame with appropriate structure
frame create results str20 theory str5 sample double mean double ci_low double ci_high double tstat double pval double dstat double pval_dist
frame results {
    set obs 16
}

est clear
    **** tests of treatment equality ****
qui forvalues i= 1/4 {
preserve
keep if treatment == `i' & experiment == 1

* Define sample label based on treatment
local sample_label = cond(`i'==1, "a", cond(`i'==2, "b", cond(`i'==3, "c", "d")))

* Get mean for roll1 (first occurrence)
regress roll1
local mean_roll1 = e(b)[1,1]

ttest roll1==3.47
estadd scalar mR1 =r(mu_1)
estadd  scalar pR1JD=abs(r(p) )
estadd  scalar tR1JD=abs(r(t) )
estadd  scalar dfR1JD=abs(r(df_t) )
local   df 		= r(df_t) 			
local critical_value = invt(`df', 0.975) 
estadd  scalar ci95hR1JD  	=r(mu_1)+ `critical_value'*  r(se) 
estadd  scalar ci95lR1JD  	=r(mu_1)- `critical_value'*  r(se)
local ci_low_roll1_jd = r(mu_1)- `critical_value'*  r(se)
local ci_high_roll1_jd = r(mu_1)+ `critical_value'*  r(se)
local tstat_roll1_jd = abs(r(t))
local pval_roll1_jd = abs(r(p))

mgof roll1 = (2/36)*roll1 + 1/36 , mc ksmirnov  reps(100000)
estadd  scalar pR1JDdis=r(p_ksmirnov) 
estadd  scalar dR1JDdis=r(ksmirnov)
local dstat_roll1_jd = r(ksmirnov)
local pval_dist_roll1_jd = r(p_ksmirnov)

* Store in results frame - first occurrence (rows 1-4)
local row1 = `i'
frame results: replace sample = "`sample_label'" in `row1'
frame results: replace mean = `mean_roll1' in `row1'
frame results: replace ci_low = `ci_low_roll1_jd' in `row1'
frame results: replace ci_high = `ci_high_roll1_jd' in `row1'
frame results: replace tstat = `tstat_roll1_jd' in `row1'
frame results: replace pval = `pval_roll1_jd' in `row1'
frame results: replace dstat = `dstat_roll1_jd' in `row1'
frame results: replace pval_dist = `pval_dist_roll1_jd' in `row1'

ttest roll1==2.5
estadd  scalar pR1FH=r(p) 
estadd  scalar tR1FH=r(t) 
estadd  scalar dfR1FH=abs(r(df_t) )
local   df 		= r(df_t) 			
local critical_value = invt(`df', 0.975) 
estadd  scalar ci95hR1FH  	=r(mu_1)+ `critical_value'*  r(se) 
estadd  scalar ci95lR1FH  	=r(mu_1)- `critical_value'*  r(se)
local ci_low_roll1_fh = r(mu_1)- `critical_value'*  r(se)
local ci_high_roll1_fh = r(mu_1)+ `critical_value'*  r(se)
local tstat_roll1_fh = r(t)
local pval_roll1_fh = r(p)

mgof roll1 = 1/6 , mc ksmirnov  reps(100000)
estadd  scalar pR1FHdis=r(p_ksmirnov) 
estadd  scalar dR1FHdis=r(ksmirnov)
local dstat_roll1_fh = r(ksmirnov)
local pval_dist_roll1_fh = r(p_ksmirnov)

* Store in results frame - second occurrence (rows 5-8)
local row2 = `i' + 4
frame results: replace sample = "`sample_label'" in `row2'
frame results: replace mean = `mean_roll1' in `row2'
frame results: replace ci_low = `ci_low_roll1_fh' in `row2'
frame results: replace ci_high = `ci_high_roll1_fh' in `row2'
frame results: replace tstat = `tstat_roll1_fh' in `row2'
frame results: replace pval = `pval_roll1_fh' in `row2'
frame results: replace dstat = `dstat_roll1_fh' in `row2'
frame results: replace pval_dist = `pval_dist_roll1_fh' in `row2'

ttest roll2==1.53
estadd scalar mR2 =r(mu_1)
local mean_roll2 = r(mu_1)
estadd  scalar pR2JD=r(p) 
estadd  scalar tR2JD=r(t) 
estadd  scalar dfR2JD=abs(r(df_t) )
local   df 		= r(df_t) 			
local critical_value = invt(`df', 0.975) 
estadd  scalar ci95hR2JD  	=r(mu_1)+ `critical_value'*  r(se) 
estadd  scalar ci95lR2JD  	=r(mu_1)- `critical_value'*  r(se)
local ci_low_roll2_jd = r(mu_1)- `critical_value'*  r(se)
local ci_high_roll2_jd = r(mu_1)+ `critical_value'*  r(se)
local tstat_roll2_jd = r(t)
local pval_roll2_jd = r(p)

mgof roll2 = (2/36)*(5-roll2) + 1/36 , mc ksmirnov  reps(100000)
estadd  scalar pR2JDdis=r(p_ksmirnov) 
estadd  scalar dR2JDdis=r(ksmirnov)
local dstat_roll2_jd = r(ksmirnov)
local pval_dist_roll2_jd = r(p_ksmirnov)

* Store in results frame - third occurrence (rows 9-12)
local row3 = `i' + 8
frame results: replace sample = "`sample_label'" in `row3'
frame results: replace mean = `mean_roll2' in `row3'
frame results: replace ci_low = `ci_low_roll2_jd' in `row3'
frame results: replace ci_high = `ci_high_roll2_jd' in `row3'
frame results: replace tstat = `tstat_roll2_jd' in `row3'
frame results: replace pval = `pval_roll2_jd' in `row3'
frame results: replace dstat = `dstat_roll2_jd' in `row3'
frame results: replace pval_dist = `pval_dist_roll2_jd' in `row3'

ttest roll2==2.5
estadd  scalar pR2FH=r(p) 
estadd  scalar tR2FH=r(t) 
estadd  scalar dfR2FH=abs(r(df_t) )
local   df 		= r(df_t) 			
local critical_value = invt(`df', 0.975) 
estadd  scalar ci95hR2FH  	=r(mu_1)+ `critical_value'*  r(se) 
estadd  scalar ci95lR2FH  	=r(mu_1)- `critical_value'*  r(se)
local ci_low_roll2_fh = r(mu_1)- `critical_value'*  r(se)
local ci_high_roll2_fh = r(mu_1)+ `critical_value'*  r(se)
local tstat_roll2_fh = r(t)
local pval_roll2_fh = r(p)

mgof roll2 = 1/6, mc ksmirnov reps(100000)
estadd  scalar pR2FHdis=r(p_ksmirnov) 
estadd  scalar dR2FHdis=r(ksmirnov)
local dstat_roll2_fh = r(ksmirnov)
local pval_dist_roll2_fh = r(p_ksmirnov)

* Store in results frame - fourth occurrence (rows 13-16)
local row4 = `i' + 12
frame results: replace sample = "`sample_label'" in `row4'
frame results: replace mean = `mean_roll2' in `row4'
frame results: replace ci_low = `ci_low_roll2_fh' in `row4'
frame results: replace ci_high = `ci_high_roll2_fh' in `row4'
frame results: replace tstat = `tstat_roll2_fh' in `row4'
frame results: replace pval = `pval_roll2_fh' in `row4'
frame results: replace dstat = `dstat_roll2_fh' in `row4'
frame results: replace pval_dist = `pval_dist_roll2_fh' in `row4'

est sto a`i'
restore
}


* Format all numeric variables to 3 decimal places and save to LaTeX
frame results {
    format mean ci_low ci_high tstat pval dstat pval_dist %9.3f
    
    * Prepare multirow for theory column
    replace theory = "JD: R_1=3.47" in 1
    replace theory = "" in 2/4
    replace theory = "FH: R_1=2.50" in 5
    replace theory = "" in 6/8
    replace theory = "JD: R_2=2.53" in 9
    replace theory = "" in 10/12
    replace theory = "FH: R_1=2.50" in 13
    replace theory = "" in 14/16
    
    * Rename variables for better column headers

    texsave using "output/table2.tex", replace ///
    frag 
}


* View the results frame
frame results: list



exit





